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Q 
U 

■ Bose-Einstein condensates (BEC) and coupled Fermi-Bose condensates near Fesh- 



Mean-field approach has recently been used to model coupled atom-molecular 



bach resonance. Sweeping of magnetic field across the resonance gives a new (non- 
linear) version of Landau-Zener problem. We investigate the structure of the corre- 
sponding classical phase space and calculate change in the action which corresponds 



(N 
> 

. to finite-rate efficiency of the sweep. We consider the case of non-zero initial action, 



which corresponds to some finite initial molecular fraction. 
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INTRODUCTION 



Adiabatic invariance is very important issue in quantum mechanics . Relation be- 
tween slow quantum transitions and change in the adiabatic invariant of a linear oscillator 
been .ud.d . g D.„^ic. of Bo.e-E..e. condensate, introduce, 
a paradigm of nonlinearity into quantum systems. In rnany mean-field models related to 



BEC physics (like nonlinear Landau-Zener models 

trapping [l^, etc), nonlinear effects that are common to classical nonlinear systems have 



macroscopic quantum self- 
classical nonlinear systems have 
assical adiabatic theory is destruction of 



13| which is encountered in different fields 



been revealed. A conceptual phenomenon of c 
the adiabatic invariance at separatrix crossings 
of physics. It is of great importance for BEC physics because change in the classical action 
of a nonlinear two-state model corresponds to probability of transition between the two 
states (modes). Here we consider a nonlinear mean- field model of a slow sweep through 
a Feshbach resonance in a quantum gas of fermionic atoms coupled to BEC of diatomic 
molecules [l^ (for brevity, we call this system a Bose- Fermi condensate). A number 
of closely related non-statio nary problems have come up recently in context of coupled 



Bose-Fermi condensates 



14 



15 



Iq ) and coupled atom- molecular BEC [171, 



mean-field approach to such problems is very interesting and not at all trivial 
the following Hamiltonian describing a fermion-boson condensate is considered: 



M). The 

1^. In [lel 



^ = E ^J^laCja + ^b^b + gY, {b^cj^cj^ + bc]/,i) , (1) 
i.o- j 

where ej are the single-particle energy levels and the operators cj^ {cja) create (anni- 
hilate) a fermion of one of the two species o" =| or | in an orbital eigenstate of energy 
Ej. In case the single-particle potential is translationally invariant, \j t) = |k t) ^ind 
|j i) = I — k I) [16]. Operators b'^ (b) create (annihilate) quanta of the bosonic field. 

The mean-field approximation of ([T]) amounts to treating the bosonic field classically, 
i.e. replacing operators 6^ and b with c-numbers in the Heisenberg equations of motion. 
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This procedure is justified provided the bosonic mode is macroscopically populated. As 
shown in [l^, in this approximation the dynamics of ([1]) coincides with that of a clas- 
sical Hamiltonian system (where classical dynamical variables are the time-dependent 
quantum-mechanical expectation values (cjjCjf), (b), and (I]o-c]^Cjo-) ). The mean-field 
approximation is also equivalent to replacing operators with classical variables and their 
commutators with Poisson brackets. It is also shown in [16] that the fermion-oscillator 
model ([T]) can be viewed as a generalization of the Dicke (Tavis-Cummings 20|) model of 
Quantum Optics: the latter model corresponds to zero fermionic bandwidth limit of (Q, 
i.e. to the case when all single-particle levels ej are degenerate, ej = /i. To demonstrate 



this, 
2kf 



16| reformulates ([T]) as a spin-oscillator model using Anderson's pseudospins (i.e., 
1, K~ = CjiCji, Kj' = cj|cjj^, where Uj = cj^Cjo-), and obtains the 



^3 

Hamiltonian 



n~l 



n-l 



j=0 j=0 

In the zero bandwidth ("degenerate") limit the Hamiltonian 
model 



(2) 

reduces to the Dicke 



Hoicke = 2/iT; + ub^b + g(b^f_ + 6f+ 



(3) 



describing an interaction of a single collective spin T = J2j with a harmonic oscillator 



Mean-field solution of (131) was discussed in 



21 



22|. The more general many-body problem 



( 1T1121) also turns out to be integrable. Explicit solutions for the mean-field dynamics of 
the model (1T1121) were constructed in [l^. Later, comprehensive solutions for the mean- 
field dynamics were obtained in [l^ using a method of separation of variables 23 



24 



25| 



which allowed to derive a complete set of integrals o 



motion 



'or (1T112|) . Quantum solutions 



of ([T]) can be obtained by the Bethe ansatz [26|, |27|. In [28|, quantum model ([T]) and a 
more general version which includes s-wave scattering interactions were solved using the 
boundary quantum inverse scattering method (QISM) as developed by Sklyanin 29| ; 
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interesting enough, through the exact solution, the spectrum can be mapped into that of 
a Schro dinger equation. 

In the present paper, we deal only with the mean-field dynamics. The mean-field 
solutions of [1^ describe dynamics of the system that has been prepared in a nonequi- 
librium state at t = 0. To model Feshbach resonance passage converting Fermi atoms to 
Bose molecules, one may use the model ([T]) with time-dependent uj. For example, in 14| . 
the model ([T]) with zero fermionic bandwidth and time-dependent molecular energy was 
considered using another approach (similar to that of [17]). To be more precise, IJ] con- 
siders the Hamiltonian H = I]k,cr ^kCk o.Ck,(j + i^(^)&o^o + 9 (Z^k Ck,tC-k,j,^o + H.c.^ , with 
Ck = being the kinetic energy of an atom with mass m, in the degenerate limit 

where = e for all k. Introducing operators J_ = ° , J- 



(Af/2)3/2 



(Ar/2)3/2 



^Ck,<T is the conserved total number of 
particles, one gets (after certain rescalings) IJ] the Heisenberg equations of motion for 
the association of a quantum-degenerate gas of fermions, 
d_ 

d 
d_ 

which depend on the single parameter 6{t) = {£{t) — 2e)/\/Ng. Mean-field limit of 



'Jx 

-Jy 
-J. 



-5{r)Jy 



I) [J: 



V2 

N 



{l + Jz 



(4) 



is obtained by replacing Jx,y,z with their expectation values u,v,w, and omitting the 
quantum-noise term ^ [1 + J^] (which is justified since the mean-field approximation is 
valid only up to terms of order 1 / \/N ; the mean-field approach to Feschbach resonance 
passage is discussed in [30, |3l|). The resulting system of equations (-it = —6{t)v, v = 
6{t)u + ^{w — l){3w +1), w = \plv) is analyzed in detail in the present paper (similar 



models arise in two-mode approximation for coupled atom-molecular BEC 



17 



18 



also, a three-mode model considered in 
portraits). 



19|; 



32] at certain conditions has very similar phase 
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We consider the model of 



ij], concentrating on the case of non-zero initial molecular 



fraction. Within the model, change in the action at the resonance passage gives the 



remnant atomic 
exponential law 



raction as a power-law of a sweeping rate parameter (instead of the 



33j of Landau-Zener linear model ). Classical adiabatic theory [1 



provides a method to calculate change in the actions in nonlinear systems. It is a well- 
known fact that action is an adiabatic invariant in a Hamiltonian system that depends on 
a slowly varying parameter. This result is based on the possibility of averaging over fast 
motion in the unperturbed [frozen at a certain parameter value] system. The situation is 
somehow different if the unperturbed system has separatrices on its phase portrait. As 
the parameter varies, the separatrices slowly evolve on the phase portrait. In particular, 
the area surrounded by a separatrix may change. Hence, a phase trajectory of the exact 
system may cross the separatrix of the frozen system. On the separatrix, the period of 
motion is equal to infinity. This results in breakdown of the averaging method and in this 
case more accurate study is necessary to describe behavior of the action. It turns out that 
at the crossing a quasi-random jump in the value of the adiabatic invariant occurs. The 
asymptotic formula for this jump in a Hamiltonian system depending on a slowly varying 
parameter was obtained in 3^, l35|, l36|] . Later, the theory of adiabatic separatrix crossings 
was also developed for slow-fast Hamiltonian systems 371], volurne preserving systems 38| . 
and was apphed to certain physical problems (see, for example, [39|, |40|, |4l|). In this paper 



we use the methods of 



361 ] to obtain a variation of the adiabatic invariant at the separatrix 



crossing. A separate paper [llj considers four different nonlinear two-mode models related 
to BEC physics, concentrating on classical non-adiabatic phenomena, and predicts new 
nonlinear effects there. Section II of the present paper introduces the model and discusses 
the structure of the classical phase space, while Section HI contains calculations of the 
change of the adiabatic invariant at the resonance passage. 
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II. MAIN EQUATIONS AND PHASE PORTRAITS 

We consider the classical Hamiltonian 



H = 6{t)w + (1 — w)y/l + w cos (j) (5) 

with a slowly varying parameter 5. Such model appears as a mean-field limit of two- 
mode approximation of coupled atomic and diatomic-molecular Bose-Einstein condensates 



Ha, 



18 



191] and in the mean-field treatments of coupled degenerate gas of Fermi atoms and 



BEC of molecules, as discussed in the Introduction. In the latter case, Feshbach resonance 
passage can be modelled by sweeping the parameter 6, and corresponding change in the 
classical action is related to the remnant atomic fraction after the sweep. 

The system can also be investigated using the equations of motion for a generalized 



Bloch vector 14l |: 



u = —S(t)v, 

V = 6{t)u + -^{w - l){3w + 1), (6) 
w = V2v, 

where the dot denotes time derivative, r = et, < e <^ 1 is a small parameter, and 
6{t) is a slowly varying parameter corresponding to the [scaled] detuning. Equations of 
motion restricted to surface are equivalent to equations of motion in Hamiltonian 
system with the Hamiltonian ([5]), where canonically conjugated variables are w and 0, 
= arctan(t>/M). 

First, we consider equations ([6]) at frozen values of parameter 6. In this case, these 
equations possess two first integrals. One of them is u'^ + v'^ — ^{w — 1)"^ {w + 1) = const. 
Corresponding to the conservation of single-pair atom-molecule coherence, the constant 
should be taken equal to zero. The equation 
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FIG. 1: Geometry of the generalized Bloch sphere (the surface + f ^ = ^{w — l)'^{w + 1)). 

+ -^{w -l)\w + l) = 0. (7) 

defines a surface of rotation around w— axis with a singular [conical] point at (0,0,1) (see 
Fig.l). Another integral of ([6]) at a frozen value of 6 is given by 

u H — j=w = const, (8) 
v2 

and at different values of the constant defines a family of planes parallel to f — axis. The 
angle between these planes and w— axis depends on the value of 6. Hence, trajectories of 
frozen system ([6]) are given by intersections of surface ([7]) and planes ([8]). At < y^, 
plane u = 5(1 — w)/^/2, passing through the singular point (0,0,1) defines a singular 
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a) 8 = -10 b) S=-2"' c) S = -l 




FIG. 2: Phase portraits of the system with Hamiltonian ([5]) frozen at different values of param- 
eter 5. 

trajectory on surface ([7]). 

Consider phase portraits of the system with Hamihonian ^ frozen at different values 
of parameter 5 [see Figure 2]. In these portraits, points (0,0,1) and (0,0,-1) of the {u,v,w) 
space are represented as segments w = 1 and w = —1 correspondingly. Point (0,0,1) is 
always a stable point of ([6]), but formally speaking this is not true for the points of the 
segment w = 1 on phase portraits of the system with Hamiltonian ([5]). If 5 is negative 
and 1^1 > v^, there is only one stable elliptic point on the phase portrait, at = and 
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w not far from —1 [see Figure 2a]. It may seem that there are also two singular points at 
w = —1,0 = vr/2, 37r/2, but they are just an artefact of the representation in variables 
w, (p; the phase trajectory passing through these points corresponds to the trajectory of 
([6]) passing through (0,0,-1), and the period of motion on this trajectory is finite. 

At 6 = —\pl a bifurcation takes place, and at > 5 > — -\/2 the phase portrait 
looks as shown in Figs. 2c, d. There are two saddle points at w = l,cos0 = 5 / \/2 and a 
newborn elliptic point at = tt. The trajectory connecting these two saddles corresponds 
to the singular trajectory on the surface ([7]). On the phase portrait of the system with 
Hamiltonian ([5]) it separates rotations and oscillating motions and we call it the separatrix 
of the frozen system. Below, we will consider " shifted" Hamiltonian E = H — 6; E equals 
zero on the separatrix. 

At 5 = the singular trajectory on the surface (JTj) is passing through (0,0,-1), and 
correspondingly on the phase portrait the segment w = —1 belongs to the separatrix 
[Figure 2 e] . At < 5 < the phase portrait looks as shown in Figs. 2 f ,g. Finally, at 
large positive values of 6, again there is only one elliptic stationary point at = vr, and 
w close to —1. 

It may be easier to understand the transition between Figure 2 h and Figure 2 i 
using also the {u,v,w) representation (see Figure 1). The elliptic point in Figure 2 h 
corresponds to the elliptic point on the surface (JTj). This latter elliptic point is close to 
the pole (0, 0, —1) of the surface ([7j). As parameter 6 grows, this elliptic point approaches 
to the pole. Consider a closed trajectory of small enough diameter on the surface ([7j), 
surrounding this elliptic point. This trajectory is given by an intersection of surface ([7j) 
and one of the planes ([H]) defined by a certain value of the constant. While S is not large 
enough, this trajectory does not embrace the pole. The corresponding phase trajectory 
in Figure 2 h surrounds the elliptic point. At 6 large enough, the trajectory on surface 
([7]) corresponding to the same value of the constant in ([8]) embraces the pole. Along the 
corresponding phase trajectory in Figure 2 i, value of (p varies from to 2tt. 
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Consider now a phase trajectory on a phase portrait frozen at a certain value of 5. 
If the trajectory is closed, the area S inside of it is connected with the action I of the 
system by a simple relation S = 2711 . If the trajectory is not closed, we define the action 
as follows. If the area S bounded by the trajectory and lines = 1, = 0, = 27r is 
smaller than 27r, we still have S = 2ttI . If 5* is larger than 27r, we put 2txI = Atx — S. 
Defined in this way, / is a continuous function of the coordinates. 

In system with Hamiltonian ^ with r = et, the action I is an adiabatic invariant of 



motion (see, for example, 13|). Far from the separatrix, it undergoes oscillations of order 
e. 

In the following, we will use the so called improved adiabatic invariant J = I + 



eu{w, 0, r) rather than /. The function u is defined as follows [see, for example, 361]]: 



u = u{w,^,t) = — 1) —dt. (9) 



27r7o V2 J dr 

The integral in (Q is computed along the unperturbed trajectory passing through the 
point [w, (f)) at the time moment t = 0; T is the period of motion on this trajectory. If 
the trajectory is not closed, T is the time necessary for a phase point on this trajectory 
to cover the distance of 27r along the 0-axis. [Definition ([9]) is valid provided that on this 
trajectory the action / and the area on the phase portrait are related simply as 5 = 27r/; 
in the following calculations of the jump in the adiabatic invariant, we need to consider 
only such trajectories.] Far from the separatrices, variation of J along a phase trajectory 
is of order e^. 

As the slow time r grows, the area bounded by the separatrix S*{t) slowly changes. On 
the other hand, a value of the adiabatic invariant associated with a certain phase trajectory 
stays well-preserved. Accordingly, phase trajectories of ([5]) can cross the separatrix. Let 
initially 5 be large in magnitude and negative, so that the phase portrait is similar to 
one shown in Figure 2a. Consider a trajectory rotating close to the singular point on 
surface ([7j). Along the corresponding phase trajectory on the plane {(p.w), value of w 
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is close to 1. We assume 1 — w on this trajectory to be small, yet finite. Hence, the 
initial value J_ of improved action J is also small. As the time goes, value of parameter 
S grows, and at 5 = —\/2 the separatrix loop appears. The area S*{t) surrounded by 
the separatrix grows with time, and the action associated with the phase trajectory stays 
approximately constant. In the so called improved adiabatic approximation this action is 
conserved at J = J_, and at the slow time moment t = such that S*{t^) = 27rJ_ the 
phase trajectory crosses the separatrix. Phenomena that take place at such crossing are 
considered in the following section. 



III. VARIATION OF THE ADIABATIC INVARIANT AT THE SEPARATRIX 

CROSSING 

First, consider the motion in the frozen system along a phase trajectory close to the 
separatrix. Let on this trajectory E = h, \h\ -C 1. [We consider such values of 6 that 
h is positive outside the separatrix and is negative inside.] The main part of the time a 
phase point on this trajectory spends in a neighborhood of the saddle points at w = 1. 
Linearizing the system near the saddle points, we find that in the main approximation 
the period of motion along this trajectory is 

T = \n\h\+b. (10) 

Here 6 is a smooth function of 6; its form is not important for the rest of the argument. 
In agreement with the formula T = 27:81/ dh we find for the adiabatic invariant I in the 
main approximation 

27iI = S*- , ^ h\n\h\ + (b + ^^=)h. (11) 

Now compute the function -u at a point of the vertex bisecting the angle between 
incoming and outgoing separatrices of the saddle point (see Figure 2c). From ([2]), ffTUl) 
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one obtains in the main approximation: 

2711^ = —==\n\h\ - b-^ + d, e = ^. (12) 

Here c? is a smooth function of 6; its form is not important for the rest of the argument. 
Consider now the separatrix crossing in the exact system with Hamiltonian E. On the 
relevant interval of slow time we have 6(t) > [the area inside the separatrix loop 
grows]. Initial values of the Hamiltonian and the improved adiabatic invariant at r = r_ 
are J_. The phase point rotates close to the separatrix. As the area bounded by 
the separatrix grows, the point comes closer to the separatrix with each turn. Denote by 
hn, In, Un values of the corresponding functions at time moments r„ when the trajectory 
crosses the vertex bisecting the angle between incoming and outgoing separatrices of the 
saddle point C outside the separatrix loop (see Figure 2c). We enumerate Tn as follows: 
To is the time of the last crossing of the vertex before crossing the separatrix, other r„ 
have negative numbers Tq > r_i > ... > r_7v ^ T-- Here ^ 1 is a large integer. Its 
exact value does not influence the result in the main approximation. After crossing the 
separatrix, the phase point continues its motion inside the separatrix loop. As the loop 
grows, the phase point goes deeper and deeper inside, rotating around the elliptic point. 
At time moments t„, 1 < n < A^, the trajectory crosses the vertex bisecting the angle 
between incoming and outgoing separatrices of the saddle point C inside the separatrix 
loop. The corresponding values of h, I, u are hn, In, Un- 
it follows from (fTO!) - (fT2|) . that in the main approximation the following expressions 
are valid: 



hn+l 
T~n+l 



2^2^ 



[ln|/i„| + 31n|/i„ — eQ\] + eb, 



(13) 
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2 2 

27r/„ = S*{to) - e(ro - r„) /^nln|/in| + {b + )h, 

1 2 - 51 J2 - 5i 



2TTUn = — 1^ \n\hn\ — b-- + d. 
2J2 - 52 4 



Here and below the values of O, b, and d are calculated at r = r*, the time of separatrix 
crossing in the adiabatic approximation, 6^, is also the value of 5 at r = r*. Summing 
the above expressions (1131) from n = —N to n = 0, we find the change of the improved 
adiabatic invariant before the separatrix crossing in the main approximation: 

2n{Jo - J-n) = -j^^ - 2« + 2^- 21n(v^/r(0)l , (14) 



^2 - 61 

where ^ = |/io/(£6)|, r(-) is the gamma function. At the separatrix crossing we obtain: 

27r( Ji - Jo) = -^S= [-K - ln(l - + 2(1 - i)Hl - + 2«] • (15) 
12 - 52 



For the change of J after the separatrix crossing we find 



2-k{Jn - Ji) 



51 



ln(l-0-2(l-01n(l-0-2e-21n(v^/r(l-0) • (16) 



Far from the separatrix, the variation of J is of order on time periods of order 1/e. 
Hence, to obtain in the main approximation the jump of J at the separatrix crossing, one 
has to sum up expression - (fT6|) . Thus, we find: 



27rAJ 



-2 , ln(2sin7rn. 

2-5i 



(17) 



Note, that this result is similar to one obtained for a symmetrical double well (see 35l. l36| ) . 
though the geometry here is different. According to 35|, l36|], the error of formula (ITTI) is 
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Formula fll7l) can be simplified. The separatrix is defined by equation E = 0. Thus, 
the area inside the separatrix loop can be calculated as 

S%5) = 2 J^^ ^[n -aTccos(^-jJ==j]dw (18) 

We are interested in derivative of S*{6) over 6. Differentiating the above integral over 
parameter 6 one obtains: 

- ^ 4^23^ (19) 



85 Js'^ y/x - 52 

We have 



where 6' = 86 /dr. Therefore, formula f[T7l) takes the following form: 

AJ= -^In(2sin7r0. (21) 

TT 

Here 6^ is the value of 6' calculated at time t = t^:. 

The value ^ strongly depends on initial conditions: a small of order e variation of 
initial conditions results generally in variation of ^ of order 1. Hence, this value can be 
considered as a random variable; its distribution should be treated as uniform on the 
segment (0,1) [see jsol, js^], see also 40] for numerically found distribution of ^ in a 
similar problem]. Formula fl2T|) is valid provided that ^ is not too close to the ends of the 
interval (0, 1): k^/e < ^ < 1 — k^/e, where A; is a positive constant, see 36||. The value 
A J in ( !2TI) should also be treated as random; we find its dispersion below. 

After the separatrix crossing, the phase point rotates around the elliptic point inside 
the separatrix loop and slowly drifts with this point to the bottom of the phase portrait. If 
6{t) is a monotonous function, the phase trajectory will never again cross the separatrix. 
Assume that 6{t) is a smooth function and 6 = 6-= const at r < r_, 5 = 5+ = const 
at r > r+. [It is assumed that 5_ < 0, 5+ > 0.] In other words, parameter 6 is slowly 
monotonically varying between two border values. Then dE/dr = at t = t±, and 
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action / coincides with J. Hence, formula fl2ip gives the variation in action /. Let the 
magnitudes of (5_ and 5+ be large enough. If the initial value oiw = is close enough to 
1, the corresponding unperturbed phase trajectory is an almost straight line (cf. Fig. 2 a). 
Hence, J_ ~ {1 — wJ). Similarly, ^ (— 1 — Thus, variation in action / corresponds 
to the remnant atomic fraction. In the adiabatic limit, value of atom-molecular imbalance 
is reversed at the passage, while the change in the action produce nonadiabatic correction 
to this result. A typical jump in the action is shown in Fig|21 

Dispersion of jumps in the action can be predicted using formula (pTI) : 

a' = 16s\6:)\-' [\n\2smnOd^=^^^ (22) 
Jo 6 

To check numerically the scaling of the jumps with e and the dispersion, we calculated 
bunches of trajectories with close initial conditions (see FigJl]). For the numerical calcu- 
lations, we used linear sweeping with 6' = 1, therefore the predicted value of dispersion is 
cr^ = (4/3)£:^. Numerically found coefficient is equal to 1.30, which is in reasonable ( 2% 
accuracy) agreement with theoretical prediction 4/3 = 1.3333. 

Consider now briefly the case when the external parameter 6 varies slow periodically 
between a positive and a negative values of large magnitudes. In this situation a phase 
trajectory crosses the separatrix once on each period of variation of 6. Each crossing can 
be characterized by two values, namely J and C,- Consider two subsequent crossings of 
the separatrix. Assume that J = Ji well before the first crossing, and the corresponding 
value of ^ is ^i. Let J2,^2 characterize the second crossing. In the main approximation, 
we have the following map: 



J2 = Jl + eF{^^) 

6 = ^i+£''GiJ2) modi. 



(23) 
(24) 
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FIG. 3: Typical jump of the adiabatic invariant (action) at separatrix crossing. 

where 

F{0 = ^^=ln(2sin7rO, 

1 

Iti 



Ti,T2 are values of the slow time r corresponding to the first and the second separatrix 
crossings, calculated in the adiabatic approximation; is a frequency of the fast motion. 
Suppose that the value gets a small variation A,^i. According to fl23|) . this leads to 
variation of the jump in the improved adiabatic invariant by a value eF'{^i)A^i, where 
F' denotes derivative of function F{C,). Thus, as follows from ( l24l) . the value ^2 obtains a 
variation A^2 = A^i + G'{J2)F'{^i)A^i, and 

^ - 1 ~ G'F'. (25) 
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0.0000 4- — • — < — . — • — I — . — . — . — I — I — I 

0.0000 0.0005 0.0010 

£ 

FIG. 4: Scaling of jumps in improved adiabatic invariant with e. For each point on the plot, 
we take a set of 100 trajectories with initial { at 6 = —10) values of w being closely distributed 
around w = 0.99 = 1 — 10^^. Final values of improved action were taken far from the separatrix 
(at 6 = 0). The line on the plot is a linear fit to data, and a perfect linear scaling can be 
seen. The slope of the line gives a ~ 1.14e, and the coefficient is in good agreement with the 
theoretical value \/4/3 « 1.15. 

This latter value can be used to describe the phase mixing in the system. If it is large, 
values ,^1 and ,^2 are statistically independent. This is the case when ,^1 is close enough to 
or to 1. Assume first that in the process of iterations of the map fl23p . (1241) all the values 
are statistically independent. Then the dynamics in the system results in diffusion 
of the adiabatic invariant and mixing in the phase space. Indeed, the process can be 
roughly modelled as a random walk in adiabatic invariant: the time interval between two 
subsequent steps is of order e~^, and the length of each step is of order e. Hence, total 
variation in J after N steps is y/Ne. After ~ e~'^ steps, which takes time of order e~^, 
the value of J (and hence the remnant part of the atomic fraction) will change by a value 
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FIG. 5: Stroboscopic shots of two trajectories without/with separatrix crossings ( a)/b) cor- 
respondingly), and with the same initial conditions. Parameters: e = 0.002, 6o = —10. In 

(a) , value of 6 was modulated as (5 = (5o(0.75 + 0.25 cos et), so no separatrix crossing occur. In 

(b) , 6 = (5ocoset, so multiple separatrix crossings result in diffusion of adiabatic invariant as 
described in the text. 



of order one. 

An example of such diffusion is shown in Fig. [5]d. The figure demonstrates the result of 
integration along one phase trajectory over a long period of time with 6 being modulated 
as (5 = (5ocos(r), So = —10. For the presentation, we plot the corresponding phase 
points on this trajectory at instances of time when 6 = 6q. For the comparison, we plot 
the results for a trajectory with the same initial conditions, but with 6 being modulating 
as (5 = (5o(0.75 + 0.25cos(r)), 6o = —10, i.e. without the separatrix crossings (Fig. [5^). 
The same stroboscopic shots produce smooth curve: the adiabatic invariant is eternally 
conserved in this system at sufficiently small e jisl . 

It is interesting to note that (for system with separatrix crossings) the variations in ^ 
may be correlated along certain trajectories. As it was shown in 4^ in a more general 
case, this results in existence of stable periodic trajectories and stability islands in the 
domain of the separatrix crossings. Total measure of these islands does not tend to zero 
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as e —>■ 0, yet it is small. Along a phase trajectory with initial conditions inside such 
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